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ABSTRACT 

Most globular clusters have half-mass radii of a few pc with no apparent correlation with 
their masses. This is different from elliptical galaxies, for which the Faber-Jackson relation 
suggests a strong positive correlation between mass and radius. Objects that are somewhat 
in between globular clusters and low-mass galaxies, such as ultra-compact dwarf galaxies, 
have a mass-radius relation consistent with the extension of the relation for bright ellipticals. 
Here we show that at an age of 10 Gyr a break in the mass-radius relation at ~ 10 6 M is 
established because objects below this mass, i.e. globular clusters, have undergone expansion 
driven by stellar evolution and hard binaries. From numerical simulations we find that the 
combined energy production of these two effects in the core comes into balance with the flux 
of energy that is conducted across the half-mass radius by relaxation. An important property 
of this 'balanced' evolution is that the cluster half-mass radius is independent of its initial 
value and is a function of the number of bound stars and the age only. It is therefore not 
possible to infer the initial mass-radius relation of globular clusters and we can only conclude 
that the present day properties are consistent with the hypothesis that all hot stellar systems 
formed with the same mass-radius relation and that globular clusters have moved away from 
this relation because of a Hubble time of stellar and dynamical evolution. 
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1 INTRODUCTION 

The half-mass radius of old globular clusters in the Milky Wa y 
depends only weakly on mass (e.g. van den Bergh et alj 1 19911) . 
If anything, a negative correlation bet ween radius and mass is 
found for the clusters in the outer halo ( I van den Bergh & Mackevl 
2004). Because this is also found for extra-galactic globular clusters 
(Jorda n et al J 1 20051 : iBarmbv et alj|2007l : iGeorgiev et alj|2009l) the 
mass-radius relation, or lack thereof, is an impo rtant aspect of the 
fundamental plan e relations of globular clusters JDiorgovskill 19951 : 
lMcLaughlij|200d) . 

Objects more massive than typical globular clusters, such 
as the recen tly disco vered ultra-compact d warf galaxies (UCDs, 
Iffilker et alj 1 19991 : iDrinkwater et alj 120031) . but also the most 
massive globular cluste rs, do exhibit a po s itive correlation be- 
tween radius and m ass JHasegan et alj|2005l : iReikuba et alj|2007l : 
iMieske et alj 120081) . Interestingly, the position of systems more 
massive than ~ 10 6 Mq in the mass-radius diagram coincides with 
the extension to low masses of the Faber & Jackson ( 1976) relation 
for bright elliptical galaxies (Has egan et al.l20 05). The mass-radius 
relation of stellar systems more massive tha n ~ 10 6 Mq h as been 
explained by the details of their formation (Murray 2009), where 



this was considered a deviation from the near constant radius of 
less massive systems. In this study we test the hypothesis that all 
hot stellar systems (globular clusters, UCDs and elliptical galaxies) 
had the same mass-radius relation initially and that the globular 
clusters (< 10 6 Mq) are the deviators because they have moved 
away from this relation because of dynamical evolution. 

Intuitively we can expect that low mass stellar systems are 
dynamically more evolved than massive systems because of their 
shorter relaxation time-scale. This evolutionary time-scale is of- 
ten expressed in terms of the half-mass properties of the system 
(Spitzer 1987) 



T rh = 0.138 



G l/2^1/2 lnA ' 



(1) 



where N is the number of stars, R^ is the half-mass radius, G 
is the gravitational constant, m is the mean stellar mass and A 
is the argument of the Coulomb logarithm and equals 0.02N < 
A < 0.11N depending on the stellar mass function in the cluster 
JGiersz & Heggiell994ar) . If we take the initial mass-radius relation 
to be of the form R^o oc Mq , then T r ho is an increasing function 
of Mq for all A > —1/3. Although the value of A is poorly con- 
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strained from observations, it is unlikely to be negative and we can, 
therefore, safely say that low mass stellar systems have shorter re- 
laxation times than massive systems immediately after formation. 

Here we consider the expansion of star clusters driven by mass 
loss due to stellar evolution and hard binaries and we present a 
description for the radius evolution including both effects, based 
on results of TV-body simulations (§|2]l. In §[3] we show that at an 
age of 10 Gyr a Faber- Jackson type initial mass-radius relation has 
been erased because of the expansion of stellar sytems with M < 
10 6 Mq . A summary and discussion is presented in §|4] 



2 EXPANSION OF STELLAR SYSTEMS 

We want to understand the evolution of the radius of a stellar sys- 
tem with a realistic stellar mass function in which the stars evolve 
and lose mass in time. This evolution is distinct from the well stud- 
ied and well understood beha viour of an equal-mass cluster (e.g. 
lHenonlll965 : Goodman 1984). Because we are mainly interested 
in the expansion we ignore the effect of a tidal cut-off. As we will 
show in §[3] the results explain the mass-radius relation of objects 
with M > 10 5 Mq, suggesting that tides are not very important in 
shaping the mass-radius relation of these objects. We first consider 
various stellar mass functions, ignoring the effect of stellar evolu- 
tion (§ 12. It , and then add the effect of stellar evolution in § 12.21 



2.1 Expansion driven by hard binaries 

The evolution o f equal-mass clusters has been studied in quit e 
some detail (e.g. lGiersz & Heggiej[l994al : iBaumgardt et al.ll2002h . 
To first order their entire evolution follows from the fact that 
gravitational systems have negative total energy, which causes 
them to always evolve away from thermal equilibrium. In the 
early evolution this results in a contraction of the core and this 
inevitably leads to the gravotherm al catastrophe, or core col- 
lapse JLvnden-Bell & EggletorJ 19801) . For an equal-mass lPlummeg 
( 1911) model the time of core-collapse is at T cc « 17T r hn (e.g. 
ISpitzeal 19871) . After core collapse the evolution is driven by bina- 
ries in the core that release energy to the rest of the cluster when 
they form and harden in 3-body interactions. This energy is con- 
ducted outwards by 2-body relaxation and in the absence of a tidal 
field this results in an expansion of the cluster as a whole, because 
escape of stars is inefficient. This increase of Rh happens on a re- 
laxation time-scale such that we can say Rh = £.Rh/T rri . If we 
integrate this relation from T cc to T, taking into account the Rh 
dependence in T rri (equation[TJ we find 



Rh 



RhO 
RhO 



1 + 

XT 

The) 



X[T - T cc ] 
Tho 

2/3 



2/3 



(2) 
(3) 



where xisa constant that relates to £ as \ = (3/2) £. In the last step 
we have used T cc ~ Tho/x as the integration boun dary, which is 
not strictly true. For an equal-mass cluster x ~ 0.14 I IHenonll 19651 : 
iHeggie & Fiull2003l) and, therefore, 1/x « 7.2 whereas for the 
Phrmmer model we have T cc /T r ho ~ 17 (see also Gier sz & Heggiei 
1994b). But equation ((3) describes the asymptotic behaviour of Rh 
for T > > Tcc and is therefore a useful approximation. It also fol- 
lows from equation l[3j that after T cc the evolution of Rh is inde- 
pendent of Rho' if we assume that TV and fh do not change in time 



then we can say Tho = T^-R^ /R h (equation [TJ and equa- 
tion (fSJ is equivalent to T rrl — xT. This means that after some 
time clusters evolve towards a mass-radius relation of the form 
i?h oc T 2 / 3 M -1 ' 3 , independent of the initial mass-radius relation. 
This is an important result and in §[3]we will show that it applies to 
globular clusters. 

T he presence of a mass fu nction speeds up the dynamical evo- 
lution ( Inagaki & Saslaw 1985), in the sense that core collapse hap- 
pens earlier dGurkan e t al. 2004 ) and the escap e rate of clusters 
in a tidal field is higher I lLee & Good man 1995). Here we estab- 
lish by means of direct TV-body simulations how the rate of ex- 
pansion, i.e. the val ue of y, depen ds on the mass function of the 
stars. We consider a lKroupal ( 120011) stellar mass function and vary 
[i = m max /m m i n , where m max and m m i n = 0.1 Mq are the max- 
imum and minimum stellar mass, respectivebjQ We consider values 
from n = 1 (= equal-mass) to 10 3 (= full mass function) in steps 
of a factor of 10. These values cover the relevant values of /j for 
real clusters. We model clusters with TV = 4096, 8192, 16384 and 
32768 particles and multiple runs are done for clusters with low TV 
and/or high /i to average out statistical fluctuations due to the low 
number of (massive) stars. The number of simulations was chosen 
to be max[l, (32768/TV)/(5 — logju)]. The initial density profile 
of all clusters is described by Plummer models in virial equilib- 
rium and during the simulation stars are taken out of the simulation 
when they reach 20i? v , where R v is the virial radius. The models 
are all scaled to the usual TV-body unit s (G = R v q = — 4Tvq = 1, 
where Eq is the total initial energy, Heggie & Mathiejl 19861) and 
we u se the kira integrator whi ch is part of the Starlab soft- 
ware ( IPortegies Zwart et al.12 001) to numerically solve the TV-body 
problem in time. At each time the values for Rh, fh and TV are 
recorded and T rrl is calculated using equation (0. For /x = 1 we 
us e 7 = 0.11 while for /x > 1 we use 7 = 0.02 as recommended 
bv lGiersz & Heggiei < fl994ah . 

In Fig.[TJwe show the (average) resulting evolution of T rrl for 
all 16 different initial conditions, specified by the number of stars 
TV and the width of the stellar mass function, [i. The increase of 
T r h after T ~ Tho/x ls dominated by expansion, because fh re- 
mains constant (no stellar evolution) and the number of bound stars 
does not change much. The dashed lines indicate different values 
of x an d it can t> e seen mal; the dependence of x on the mass 
function can to first order be approximated by x ~ 0.1/i 1 ' 2 (or 
X ~ 0.1(m max /m) 0,7 ). This scaling roughly recovers Henon's re- 
sult for equal-mass models. In summary, we see from Fig. Q] that 
T/T r h increases until T ~ T r ho/x> after which T/T r h ~ constant. 

In § 12.21 we repeat the simulations with (j, — 10 3 and turn 
stellar evolution on such that fi naturally decreases from 10 3 at 
T = to \i ~ 10 at T « 10 Gyr during the simulation because of 
stellar evolution. 



2.2 The combined effect of stellar evolution and binaries 

The time-scales of stellar evolution are set by the stellar interi- 
ors and are independent of the relaxation time-scale of the cluster 
wherein the stars evolve. We thus expect the details of the evolu- 
tion to depend on a combination of the stellar evolution time-scale 
and the relaxation time of the cluster. Here we show that the re- 



1 We use the ratio m max /rre m j n because it is easy to relate to real clusters. 
iGiirkan et al. (20040 show that the relevant parameter is m max /?fi which 
captures variations in m max and the slope of the mass function. 
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Figure 1. Evolution of the half-mass relaxation time, T r h, for clusters with 
different N and different fi. The ./V-body unit of time, T^yn, can be related 
to physical units through T dyn = (GM/Ft^)- 1 / 2 . A Kroupa (2001) mass 
function is used for the stars in the range 0.1 ^ iti/Mq ^ O.lp. Clusters 
of different N and the same fi evolve to the same T rh w \T after core- 
collapse (equation [3). This asymptotic behaviour is roughly matched by 
the relation \ ~ 0.1/j 1 ' 2 , shown as dashed lines. The T r h values of the 
equal mass clusters are calculated using a slightly different argument in the 
Coulomb logarithm (A = 0.11A r , equation[TJ as compared to the multi- 
mass clusters (A = 0.02N). For clarity the /i = 1 curves are only plotted 
forT>3T rh0 . 



suiting expansion still depends in a simple way on the dynamical 
properties of the cluster as a whole. 

We want to consider a large range of T r ho with our simula- 
tions to cover a parameter space that is relevant for real globular 
clusters. Because computing times limit us to N < 10° with di- 
rect A^-body simulations, we vary both N and the initial half-mass 
density, ph = 3M/(87T.R^). We consider 15 different values of 
T rh0 ranging from T rh0 « lMyr ([JV,logp h ] = [8192,6]) to 
T rh0 ~4Gyr([Af,logp h ] = [131072, 1]), with p h in M pc~ 3 . 
Here T r ho is increased by increasing A^ by factors of 2 and by 
decreasing ph by factors of 10. We again use the kira integra- 
tor and the stellar evolution package SeBa for solar metallicity 
( IPortegies Zwart et al. 2001). We use a Kroupa (2001) initial mass 
function between 0.1 M and 100 M©, which has m « 0.64 M . 
The retention fraction of black holes and neutron stars was set to 
zero. 

In Fig. [2] we show the resulting expansion in the form of 
Rh/Rho as a function of T r ho at different ages. The asymptotic 
behaviour of these runs can easily be understood by considering 
the extremes. Clusters that are dynamically young (low T/T T ho) 
expand adiabatically in order to retain virial equili brium after stel - 
lar mass loss. The continuous loss of mass from a lKroupal (2001) 
mass function together with the stellar evolution prescr iption of 
Starlab (Appendix B2 of IPortegies Zwart et all 1200 ll) leads to 
a reduction of the total cluster mass 



M w Mo 



T > T», S « 0.07, T, « 2 Myr. (4) 



T» 



(5) 



This adiabatic expansion is slow in time and gives a maximum in- 
crease of Rh/Rho ~ 2 after a Hubble time. At the other extreme 
we have clusters that are dynamically old (high T/T r ho) and they 
expand quickly in a way that is comparable to what we have seen 
in § 12.11 We propose a function that stitches together these two ex- 
tremes in an attempt to match Rh/Rho f° r all values of T/T r ho 
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(6) 



Here \t is a parameter comparable to \ of § 12.11 but now time- 
dependent due to the variation of the mass function, and its value 
at an age T is found from a fit of equation ® to the results of the 
A^-body runs. In Fig. [2] we show the fit results as full lines and the 
resulting values of \t are indicated. It shows that equation ((6) pro- 
vides a good description of the evolution of R^/R^o- The relation 
between \t and T is well approximated by a simple power-law 
function 



XT 



T* < T < 20 Gyr. 



(7) 



If we now define T» = min([2Myr, T]) we have a continuous 
function for R h (T lh0 ,T), or R h (M ,R h0 ,T) for all T. For high 
T/Trho we find from equation ® that T r h = (mo/m) 1 ' XtT oc 
T 0,74 (equations|4]&|7]l. We indicate below how the small deviation 
from a linear scaling with T (as found in § I2.lt can be interpreted 
in terms of the evolution of the mass function. In § 12. II we found 
T r h oc /j, 1 ' 2 T, If we appr oximate the ma in sequence life-time of 
stars by i ms ex m~ 25 (Br essan et alj 1993b and thus fj, oc T" 1/2,5 , 
then from the changing mass function we expect T r h oc T°' s , very 
close to what we find from the numerical simulations. We conclude 
that the evolution of the cluster is 'balanced' when the second term 
on the right-hand side of equation ® dominates (high T/T r ho), in 
the sense that the energy flux at the half-mass boundary that drives 
the expansion is provided by the production of energy in the core 
by binaries and stellar evolution combined. Stellar evolution in fact 
slows down the expansion rate, because R^ is determined by the in- 
stantaneous width of the stellar mass function, resulting in a smaller 



In this regime the radius thus evolves as (e.g. Hills 1980) 



7?h at old ages than if p had stayed constant at p = 10 (§ 12. It . 
When the first term on the right-hand side of equation ® domi- 
nates (low T/T r ho) the evolution is unbalanced and we have the 
usual adiabatic expansion . 

The fact that the interplay between dynamical evolution and 
stellar evolution is in fact quite simple can be understood from the 
energy budget. The total energy of a stellar systems depends on 
M and R^ as E oc —M 2 /R^. Together with equations lf4]l&l[5} 
we find that E evolves in time as E oc — (T/T t )~ 3 because of 
mass-loss and (adiabatic) expansion. The rate of energy change as 
a result of stellar evolution is then Ssev oc \E\/T, where the con- 
stant of proportionality depends on the degree of mass segregation: 
it will be higher when stars lose mass from the centre and/or when 
the density profile is centrally concentrated. The rate of energy in- 
crease due to binaries and relaxation is similar. This is because 
E oc \E\/T T h and T r h oc T (ii l2.U . After some dynamical re- 
laxation mass-loss by stellar evolution will be predominantly from 
the core because of mass segregation and this will boost Esev- We 
tentatively pose the idea that Esev acts as a central energy source 
that is subject to a feedback mechanism comparable to what hap- 
pens with binaries: if Esev is too high, the core expands and the 
central potential decreases and Esev drops. If Esev is too low, the 
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Figure 2. Expansion from the iV-body runs including the effect of stellar 
evolution together with the functional fits (equation[6] full lines) 



core contracts, thereby increasing the depth of the cent r al po tential 
and increasing -Esev- This fits in the view of lHenonl i 19751) that 
"the rate of flow of energy is controlled by the system as a whole, 
not by the singularity". One of the consequences is that there is no 
sharp transition between a stellar evolution dominated phase and a 
relaxation dominated phase. 



3 APPLICATION TO OLD STELLAR SYSTEMS 

With the expression for the evolution of the radius as a function 
of T r ho at hand we can easily calculate the evolution of R^ for 
any initial mass-radius relation. We apply our result to the mass- 
radius relation of old and hot stellar systems in the mass range ~ 
10 4 -10 8 M Q . 

The original Faber-Jackson relation relates the central veloc- 
ity dispersion of (bright) elliptical galaxies to their total luminos- 
ity. lHasegan et alj J2005h have converted this result into relations 
between M, R^ and surface density. The resulting mass-radius re- 
lation (their equation 15) with an additional log 4/3 to correct for 
projection is log(_R h /pc) = -3.142 + 0.615 log(M/M Q ). They 
show that this relation matches the objects with M > 10 6 Mo 
(UCDs, massive globular clusters and their dwarf-globular tran- 
sition objects, DGTOs) in the mass-radius diagram. Because this 
concerns collision-less systems, we can safely assume that 2-body 
relaxation has not affected this relation and it should, therefore, 
reflect the initial relation. To get an expression for the initial mass- 
radius relation we only need to correct for mass-loss by stellar evo- 
lution and the subsequent adiabatic expansion. For T = 10 Gyr, 
we find M/M = Rho/Rh ~ 0.55 (equations 0]&[5} and thus 



log 



RhO 

pc 



= -3.560 + 0.615 log 



M \ 
Me) 



(8) 



In Fig. [3] we show how this initial mass-radius relation evolves us- 
ing our result from equation (|6j together with data points that cover 
the mass regime we are interested in. For high T/T r ho the radius is 
set by Mo, independent of Rho, while for low T/T r ho we are seeing 
roughly the initial mass-radius relation. By construction the right- 



hand side of the 10 Gyr line coincides with the re presentation of 
the Faber-Jackson relation of lHasegan et alj ( 120051) . From solving 
dRh/dM = in equation {6]l we find that at an age of 10 Gyr the 
break between the two regimes occurs at M u fa 1.1 x 10 6 Mq and 
at that age systems with this mass have T T h/T ~ 0.8. lMieske et alj 
(2008) noticed already that the break occurs at systems with T r h 
roughly equal to a Hubble time. In this paper we give a quantitative 
explanation for it. 



4 SUMMARY AND DISCUSSION 

In this study we provide the arguments that explain why there 
is a break in the mass-radius relation of hot stellar systems at 
~ 10 6 Mq . We show that the mass-radius relation of the massive 
systems (> 10 6 Mq) is only slightly affected by stellar evolution 
and represents, therefore, approximately the initial mass-radius re- 
lation. The origin of this relation needs to be searched fo r in the 
detail s of their formation and is not discussed here (see e.g. Murray 
2009). Combining scaling relations for the (adiabatic) expansion 
of clusters because of stellar evolution with relations for expan- 
sion due to 2-body relaxation we present a simple formula for the 
radius evolution as a function of initial mass, radius and time. Ap- 
plying this result to a Faber-Jackson type initial mass-radius rela- 
tion (the representat ion in units of mass and radius are taken from 
lHasegan et alj|2005l) we show that at an age of 10 Gyr a break oc- 
curs at ~ 10 6 M©. This break can be thought of as the boundary 
between collisional systems (T r h < age) and collision-less systems 
(T rh >age). 

For young massive clusters there is also no obvious correla- 
tion between radius and mass/luminosity (Zepf et al. 1999; Larsen 
2004; Scheepmaker et al. 2007; Portegies Zwart et al. 2010). From 
Fig. [Hit can be seen that for clusters with an age of ~ 10 — 100 Myr 
there has already been significant expansion of clusters with masses 
< 10 5 Mq. Although this break mass depends on the initial mass- 
radius relation, it at least qualitatively shows that at young ages 
most clustery are affected by the expansion we consider here. It 
is worthwhile to compare the theory to the parameters of young, 
well resolved star clusters (e.g. lMackev & Gilmorell2003l) . We em- 
phasise that the balanced evolution provides a lower limit to cluster 
radii. If clusters form above the relation marked initial in Fig.[3]ihen 
they expand only slightly because of stellar evolution at young ages 
until T/Ti-ho is high enough for the balanced evolution/expansion 
to start. The mass-radius relation of young clusters is important 
in the evolution of cluster populations. This is because in the 
early evolution clusters suffer from encounters with the molecu- 
lar gas clouds from which they form. The time-scale of disrup- 
tion due to su ch encounters scales with the density of the cluster 
( ISpitzeilll958r) . If all cluster have the same density, their disrup- 
tion time-scale is independent of their mass. For a constant radius 
the time-scale of disruption becomes strongly mass-dependent be- 
cause then ph oc M and for a constant T r h we have ph oc M 2 . 
The mass-radius relation, therefore, determines the properties of 
the clusters that survive continuous encounters with massive clouds 
JGieles et alj2006l ; lElmegreenl201Cjh . 

We have ignored the tidal limitation due to the host galaxy. 
Once the density of a cluster drops below a critical value, de- 
pending on the tidal field strength, our result will overestimate 

2 The luminosity function of young clusters is a power-law distribution 
with index ~ — 2 such that a typical young cluster population only has a 
small fraction of its clusters in the massive (> 10 5 Mq) tail. 
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Figure 3. Mass-radius values for hot stellar systems. The values for 
globular clusters in the Milky Way, the Magellanic Clouds and Fornax 
are taken from McLaughlin & van der Marel 12005). The clusters in M31 
are from iDubafh & Grill maii 11997). The values for globular clusters in 
N GC 5128 (Cen A), UC Ds and DGTOs are from the compilation presented 
in iMieske et al. 12008). The lines show the evolution of the mass-radius 
relation using the Faber-Jackson relation, corrected for stellar evolution, as 
initial conditions. The break at ~ 10 6 Mq at T rj 10 Gyr is because lower 
mass objects have expanded. 



the radius of such clusters because the presence of a tidal limi- 
tation will prevent further growth. The good agreement between 
the simple model presented here and the data points suggest that 
at least to first order the positions in the mass-radius diagram of 
objects with M > fewx 10 4 M© is not much affected by tides. In- 
cluding a tidal field would bend down the curves at low masses. 
The transition from expansion dominated evolution to Roche-lobe 
filling evolution is considered in more detail in a follow-up study 
JGieles. Heggie & Zhadl20ld) . 
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